function fitness_history = pso(n_particles, n_variables, max_iter, lb, ub, objective_function)
    w_inertia = 0.6;  % 默认惯性权重
    c1 = 0.5;         % 默认个体学习因子
    c2 = 0.5;         % 默认社会学习因子
    
    % 粒子位置初始化
    X = lb + (ub - lb) * rand(n_particles, n_variables);  % 粒子位置
    V = rand(n_particles, n_variables);                   % 粒子速度

    % 初始化个体最优解和全局最优解
    Pbest = X;  % 个体最优解
    Pbest_fitness = arrayfun(@(i) objective_function(X(i, :)), 1:n_particles);  % 个体最优适应度
    [gbest_fitness, gbest_idx] = min(Pbest_fitness);  % 全局最优适应度
    Gbest = X(gbest_idx, :);  % 全局最优解

    % 记录适应度函数的变化
    fitness_history = zeros(max_iter, 1);  % 存储每一代的全局最优适应度值

    % PSO主循环
    for iter = 1:max_iter
        for i = 1:n_particles
            % 更新粒子速度和位置
            V(i, :) = w_inertia * V(i, :) + c1 * rand() * (Pbest(i, :) - X(i, :)) + c2 * rand() * (Gbest - X(i, :));
            X(i, :) = X(i, :) + V(i, :);
            
            % 限制粒子位置在搜索空间内
            X(i, :) = max(min(X(i, :), ub), lb);
            
            % 计算新位置的适应度（目标函数值）
            new_fitness = objective_function(X(i, :));

            % 更新个体最优解
            if new_fitness < Pbest_fitness(i)
                Pbest(i, :) = X(i, :);
                Pbest_fitness(i) = new_fitness;
            end
        end
        
        % 更新全局最优解
        [min_fitness, min_idx] = min(Pbest_fitness);
        if min_fitness < gbest_fitness
            gbest_fitness = min_fitness;
            Gbest = Pbest(min_idx, :);
        end
        
        % 记录每次迭代的全局最优适应度值
        fitness_history(iter) = gbest_fitness;
    end
end